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Abstract 

Harmonic grid generation methods for multiply connected plane regions and 
regions on curved surfaces are discussed. In particular, using a general 
formulation on an analytic Riemannian manifold, it is proved that these 
mappings are globally one-to-one and onto. 
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INTRODUCTION 


One of the most Important tasks of computational physics is the generation 
of boundary conforming grids in the computational domain. This problem can be 
framed mathematically as the construction of certain coordinate charts on a 
given Riemannian manifold. In this paper the method of harmonic mappings for 
two-dimensional manifolds is considered. The method of two-dimensional 
harmonic mappings has been used, for example, to generate grids in multiply 
connected plane regions [1] and to generate coordinate surfaces for three- 
dimensional domains using the Gauss equations [2] . We will see that the 
harmonic mapping method can also be used to generate grids on a portion of a 
given analytic curved surface. We will describe and analyze these methods by 

unifying them to a general mapping problem on a simply connected Riemannian 
manifold. 


2. HARMONIC MAPPING METHOD FOR AN ANALYTIC RIEMANNIAN MANIFOLD 

The simplest grid generation problem of this kind can be formulated as 
follows. Let Q be a simply connected region on a two-dimensional manifold 
M 2 with compact closure (see Figure 1). Suppose (0,40 is a given (natural) 
coordinate system on this manifold. The problem is to find a boundary 
conforming coordinate system (S,T)). This is accomplished by first mapping 
9 &q in a continuous 1-1 manner to the boundary of the rectangle Qj in the 
(S,U) plane. To map the interior of Q Q to the interior of Q we use the 
following conditions : 



( 2 . 1 ) 


AS = 0, 


An = 0 


where A is the Laplace-Beltrami operator. 

We will see later that more complicated grid generation problems can be 
treated in this manner. Let us investigate more closely this technique. We 
set (x 1 , k 2 ) = (0,t) and (x 1 ,* 2 ) = (5 ,n ) for ease of notation. Then g^ 
and £ are respectively the metric tensors of these coordinate systems. 
Then from the elementary theory of Riemannian geometry we have the following 
explicit expression for the Laplacian (Beltrami) of any scalar function <f> on 


( 2 . 2 ) 


A1 aBf-9 <f> 3<J> 1 

A+ = g L— ; 5- r *R 


3 x a 3 x^ 


3x XJ 


where a,B,X = 1,2. In this context we have used the Christoffel symbol of 
the second kind in the x a coordinates. It is well known that this symbol 
obeys the transformation laws, 

* c t 


(2.3) 


X _ 3x X 3 2 x Y f pM 3^ 3^ 9x X 

3x Y 3x a 3x P 3x a 3x S 3^ ’ 


where, of course, is the Christoffel symbol of the second kind for the 

coordinate system x • 


2*2. Formulation on a Given Manifold 

9 9 — ( 

We note here that if the manifold M = E (Euclidean space) then x 
can be chosen as the Cartesian coordinate system yielding f^. =0 [5 


- 2 - 


This is not possible for a general manifold. However, since any two- 
dimensional Riemannian manifold is conformally flat [5], it is always possible 

Q 

to choose (locally) a coordinate system x known as the isothermal 
coordinate system which satisfies 

(2.4) g ae = A* 1 ,* 2 ) V 


In this case, one can verify that 


(2.5) 


C - “ 


yielding 


A . a 3 3 <|> 

A<t> = g 


3x a 3x^ 


3 x^ 3 <f) a 3 3 ^ x^ 

i-Y 77 s 77 7.fi 


3x 9 x 


9 x 9 x 


Xf we choose <|> = x (one of the coordinate functions) then 
considerable simplification takes place and we have 


6 2 — Y 

A 6 3x° ct3 3 Z x Y 

Ax g 5 - • 

3 x 3 x 3 x 

By imposing the condition Ax° = 0, we in effect impose conditions on the 
resulting metric tensor. Yielding the equation 


ag 9 2 x X 
5 9x° 9x B 


= 0 , 


X = 1,2 


as the generating equations to be solved in the (x^x 2 ) plane. For example, 
we consider a portion of a spherical surface. In this case we can obtain an 
isothermal natural coordinate system x* using the stereographic projector 
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"“1 2 

x + ix = tan(^/2) exp(i0) 

where ( 0,10 are the spherical coordinates. Imposing Ax 5 = 0, 5=1,2 yields 
as before 


ag 3 2 x X 
3x a 3x e 


= o, 


X = 1,2 


with Dirichlet boundary conditions. A practical grid system for a conical 
wing shape [3] is given in Figure 2. Similar computations could be performed 
for other simple surfaces where some isothermal coordinates are explicitly 

known. For certain simple surfaces, isothermal coordinates are given in 
[4,5]. 


2.3. Surface Construction Using Harmonic Mappings 

We now consider the problem of generating a coordinate surface in E 3 
connecting two one-dimensional boundaries (lying in two boundary surfaces). 
Again, the idea is to place implicitly restrictions on the metric tensor (see 
Figure 3). Let the unknown manifold M 2 be immersed in E 3 with Cartesian 
coordinate system x^, i = 1,2,3. The tangent vectors are 


3x 


35 ' 


a * 


i = 1,2,3, a = 1,2 


1 2 

where (£ ) are some coordinate system for M 2 . The Gauss [7] equations 


are 




N 


i 
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where "H" denotes covariant differentiation, Q a p is the second fundamental 

Ct ft 

form and N* is the normal to the manifold. Multiplying by g yields 


a 2 x 1 _ap ^ r i B h ap . 

8 »„a :j - g r aB B \ + F hk B a B p g (k l + k 2 )N * 


ar 


P 


where and are respectively the Christoffel symbols of the second 

kind of the l a coordinate system on the manifold and the x 1 coordinate 
system on E * Also (k^ + k2> is twice the mean curvature. If x are 
Cartesian coordinates then = 0 and recalling our earlier development, we 

obtain 

a(3 5^ «i , , v, T i 

g — p + B. = (k. + k ? )N • 

5C a X 12 

If we impose A£^ = 0 for X = 1,2, we get 


ap 5^ x* /t , , NXT i 
g — 3* = (k, + k )N 

bi a ac p 1 2 


o 

which could be used as generating equations for M • This technique has been 
exploited by Warsi [2] , [6] . 


3. MULTIPLY CONNECTED DOMAINS AND THE CONSTRUCTION OF AN ABSTRACT RIEMANN 
SURFACE 

We will first describe briefly the construction of an abstract Riemann 
surface for an algebraic function as a compact and orientable surface that 
could be realized in the Euclidean space E 3 . Consider the function 
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2 2N 

w = IT (z - r. ). 

i=l 1 

Then w will be single-valued in the two-sheeted cut plane as shown in Figure 
4. Projecting each of these planes stereographically to their Riemann spheres 
and joining the spheres along the cuts after appropriate rotations, we get a 
closed surface that is topologically equivalent to a surface of genus N-l. 
We have thus obtained a surface on which w is single-valued, in a form which 
could be realized in E^. 

If we are now interested in generating grids by harmonic mapping on a 
mulitply connected domain with appropriate cuts as shown in Figure 5, we 

see that this domain will have a corresponding image on a suitable abstract 
Riemann surface S. We will now introduce the Riemannian metric in this 
surface to convert it to a metric Riemann surface (conversely, an analytic 
Riemannian manifold with isothermal coordinates will become an abstract 
Riemann surface). Thus, we see that the problem of generating grids on a 
multiply connected plane domain and that on a portion of an analytic curved 

A 

surface are equivalent. Let S be the universal covering surface (simply 
connected) of our Riemann surface S. On this covering manifold the domain of 
interest will have many (disconnected) images and these images can be 

transformed to each other by what is called the covering transformation. We 
will thus work with the fundamental domain on this manifold that is left 
invariant by the group of covering transformations. Since any simply 
connected Riemann surface can be conformally mapped to a plane [4], the 

A 

covering surface S will have a one-to-one conformal transformation with a 
domain G on the plane (see Figure 6) . 


Let P be a point on the Riemann surface S and P be one of its images 

A 

on the universal cover S. The projection mapping ^ is defined by 

7i(P) = P, 

i 

(note it is not single-valued)* Let f map S to G, then f is 

univalent and 

—1 —2 

a = x -fix = f (P) 
or 

p = (n o f -1 ) (a) =^(a) . 

Thus, cr = ^~^(P) is the uniformizer for the Riemann surface S. This 
procedure has provided us with a global family of isothermal coordinates 
(x ,x ) on the manifold S. We have to note that x is not single-valued 
if S is not simply connected and the potential o may have stationary 

points (for a surface of genus g, a potential with m poles will have 

2m + 2g - 2 stationary points). However, it is to be expected that a 

suitable o can be found with no singularities in the domain of interest. 

We will now consider the problem of producing a grid on Q as above using 
harmonic mapping techniques. This method will provide non-isothermal 
coordinates in general. However, we will use the fact that the domain of 
interest Q on the manifold could be covered by a global isothermal 
coordinate system to develop the rest of the theory. Thus, the harmonic 
mapping problem on a curved surface or a multiply connected plane domain with 
multivalued coordinates (x^,x^) has been reduced to a problem of finding 
single-valued coordinates (x (x ,x ), x (x ,x )J on a simply connected 
domain in the uniformizing plane G. 
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We will now express the Laplace-Beltrami operator in these uniformizing 


variables, 


- 4 ( 


2 v. _i 2 

SCxV 




), 


X - X (x^ ,x^) , 


Thus, we have Ax a = 0, a = 1,2, in the image Qq on the plane region G. 

Let us now show that this grid generation procedure does map the interior 
of to the interior of in a one-to-one manner. The problem has been 
reduced to the uniformizing plane and therefore the following theorem applies 
[8], we include the proof for completeness. 


THEOREM 1: The coordinate functions (x^,x^) obtained by harmonic 

mapping techniques map into the rectangle have non-vanishing 

gradients in 

Proof : Suppose Vx 1 =0 at a e then 

1 1 * 

W (a) = x + i x , 


where 


, 1 * 


is an harmonic conjugate of 


X 1 , 


must satisfy 


aw 


0 


at 


which implies that Z (a) = W(a) - W(cFq) should have a zero of order greater 
than or equal two. If this is true, then the argument of Z(a) around afl^ 

is at least 4ir which means that x*(a) - x*(<Jq) should vanish at least four 
times on 8^. But we chose to be a rectangle and therefore this 
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function has exactly two zeroes on 3 ^q. Therefore, we conclude that Vx* 

2 ° 

and Vx are never zero in We remark that if had been a more 

complicated region, then this proof would be no longer valid. 


THEOREM 2: Let T be the harmonic mapping of into That is, 


T : 3^0 -»• 3J2j homeomorphically 

and 

o 

AT = 0 on 0 . 

— 0 


Then if T _is C 1 _on figMPiJjLi where {Pj}^ C 3fl ? , then the Jacobian 

o 

of T never vanishes in 0 . 

Q 


Proof: Let 


p e a^ 0 \{p i > 


N 

i=l * 


Then 30 q is smooth at P and we may unambiguously define an orthogonal 
coordinate system s(P), n(P) where F(P) is the tangent vector on 3ftg 
obtained by traversing 351^ at constant velocity and rT(P) points into 0^, 
This yields an orientation of 38 q. Now let ^ R n ^=l ^ 0 with -»■ p. 

Then computing the Jacobian of T at R n in the coordinates s’(P), ri(P), we 
have for C * 0, 

C.(x* x^ - x^ x*) (R ) = J(T)(R ). 
s n s n n n 

Since T is C*, these values converge to the value at P, namely 
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C.fx^ x 2 - x 2 x*)(P). 

11 s n s n' 

But since P is a boundary point and is a rectangle, then one of the two 

terms vanishes and the remaining term is nonnegative if both curves have been 
oriented counter-clockwise. Thus, we conclude that the Jacobian of this 
mapping is well-defined at the boundary and nonnegative. 

Now consider 


, N (3x 2 )/Ox 1 ) + i[(3x 2 *)/(3x 1 )] 

x(a) = 1 =3 1* =TT ’ 

(3x )/(3x ) + i( (3x )/(3x )) 

o ^ o 

where x is the conjugate harmonic function to x • Thus, making use of 
the Cauchy-Riemann equations we have 

Re(x<°)) l^‘| 2 , 

3 (x 1 ,x Z ) 

o 

which is harmonic in SJq and hence, is positive in by the maximum 

principle. Here we have used the fact that Vx 1 * 0 and the nonnegativity of 
the Jacobian in the boundary. 

Thus, the harmonic mapping method produces for multiply connected plane 

domains as well as portions of analytic Riemann manifolds, a locally 1-1 

mapping function. We will now prove that these mappings are globally 1-1. 

o 

Although our present focus is on mappings in R , the next two lemmas are 
presented in R^ for future extensions to mapping problems in higher 
dimensions. 


- 10 - 


LEMMA 1: 


Let ft„ and ft, be compact subsets of 1^ with ft = ft 
0 1 0 0 

and ft ^ convex. Let T be a continuous mapping from ft q into T$ which 


o 


is differentiable in ft and satisfies 


1. T(9ft Q ) = 3ft x 


2 . 


T'(x)| * 0 


for all x e ft 

- 1 


then T maps ft^ onto ft^. 


Proof : The first thing to notice is that T is an open map on ft and 

o 

hence, T(ft ) is open. We claim that T(ft ) C ft,. If not by compactness, 
0 U l 

* ° 

there is an x eft satisfying 

- 0 


& 

0 < dist(T(x ), ftjJ = max dist(T(x), ftj) 


But T is open on ft and ft is convex, which implies that 

0 1 

y dist(y, ft^) has no local maxima. This contradicts the assumption that 

. o o 

T(ftg) (f ftj. Thus, we have T(ft ) C ft and even T(ft^) C ft^. If we fail 

o 

to have equality, there is a ye ft^\T(ft^). But for any such y, we must have 


inf IIT(x) - yll = inf IIT(x) - yll , 


* e % 


x e 9£V 


o o o 

since T is open on SI • However, one can always produce a y e SI \T(ft Q ) 

o 

arbitrarily close to T(H^) an< * relatively far from 9 Sl^ = T(9Aq) to violate 
this condition. This completes the proof of this lemma. 
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LEMMA 2: Let ftp and ft^ be compact homotopic subsets of with 

ft j convex . Let R^ , 0 < a < 1 be the homotopy with 


vv 


Rj(x) = x 


for x e 


|R'(x)| * 0 
a — 


for x e ft 


Let T^, 0 < a < 1 be the mapping from ft^ I?* satisfying, 


T^Cx) = 21 for x s ft^ 


|T'(x)| * 0 
a 


for x e ft 


T on ) = 3ft, 
a a 1 


T_ ° R_ is a homotopy on ft, . Then T is a one-to-one and onto 
a a 1 a — 

transformation from ft^ into ftp. In particular, Tq maps ftp one-to-one 
and onto ft, . 


Proof : Recall that for a differentiable map T : ft ->■ I? (for ft compact 

subset of we can define 


- 12 - 


(a) 


d = degree (T,y,fi) = £ sign|T'(x)| 

T(x) = y 
o 

xeO 


, ° 

provided y^ T(90) and for each xe 0 above, |t'(x)| * 0. Furthermore, 
this integer d is invariant under homotopy (homotopy invariance theorem) 
provided no solutions are introduced on the boundary [9 ] . Now to prove the 
lemma, set 


S = T o R . 
a a a 


Thus is a homotopy on 0^ and for each x e 0 


|S'(x)| * 0 


and hence, |S'(x)| is always positive or always negative. But |S'(x)| = 1 

and hence, |S'(x)| > 0 for all x e 0 • We know from Lemma 1 that T 
a — — 1 a 

o o 

maps 0 onto 0 and 0 onto 0 . Thus, S maps 30. onto 30, and 

o a o 1 o a a 1 1 

onto 0 r Let y e 0^ and note 

degree(S 1 , y, flj) » 1 

since is the identity. Both is a homotopy which introduces no 

solutions on the boundary and hence, by degree theory 


degree(S Q , y, flj) = 1, 

. ^ | o 

since |S (x) | > 0 for all x c J! . Thus, (a) actually counts the number of 
a - - 1 
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x's, so that 


s a (x) = y* 

O 

Since this number is one we conclude that, S maps ft , one-to-one and onto 

a 1 

o 

. and hence T maps ft one-to-one and onto ft, and finally ft n one-to- 
l * a r a 1 u 

one onto ft^ as was to be shown, 

CONCLUSION 

Grid generation methods using harmonic mappings have been analyzed in a 
unified framework and a rigorous justification for these methods is given. We 
have proved that for these mapping methods, mesh intersection, or overspill 
cannot occur. Our conclusions are limited to two-dimensional manifolds, 
although substantial portion of the theory developed is applicable for higher 
dimensions. An extension of this theory should include the generalization to 
higher dimensions and to the method of nonharmonic mappings. 
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